/** blocksorting.d - A implemention of Burrows-Wheeler transform

  Written in the D programming language 1.0

Authors:	Wei Li <oldrev@gmail.com>
Copyright:	Copyright (c) 2007 by Wei Li.
License:	BSD
*/

module dotmars.compression.algorithm.blocksorting;

// begining of MTF coding
////////////////////////////////////////////////////////////////////////////////

private template MtfCodingBase()
{
	const int SymbolMax = 255;

	private ubyte[SymbolMax + 1] m_order = void;

	void init()
	{
		for(size_t i = 0; i < m_order.length; i++)
			m_order[i] = i;
	}
}


struct MtfEncoding
{
	mixin MtfCodingBase;

	ubyte encode(ubyte c)
	{
		int i = 0;
		for(; i < m_order.length; i++)
		{
			if(m_order[i] == c)
				break;
		}

		for(int j = i; j > 0; j--)
			m_order[j] = m_order[j - 1];

		m_order[0] = c;

		return cast(ubyte) i;
	}
}


struct MtfDecoding
{
	mixin MtfCodingBase;

	ubyte decode(ubyte i)
	{
		ubyte c = m_order[i];

		for(int j = i; j > 0; j--)
			m_order[j] = m_order[j - 1];

		m_order[0] = c;
		return c;
	}
}

////////////////////////////////////////////////////////////////////////////////

// BWT core 
// The original C source of "bwtsort.c" has written by Karl Marlbrian
// Reference:	http://www.wxwidgets.org/wiki/index.php/Development:_Karl_Marlbrian_BWT
struct BlockSorting
{
	//  these two values are stored together
	//  to improve processor cache hits

	struct KeyPrefix
	{
		size_t prefix;
		size_t offset;
	}

	//  offset/key prefix
	//  for qsort to use

	private KeyPrefix[]	m_keys;
	private size_t[]	m_rank;

	//  During the first round which qsorts the prefix into
	//  order, a groups of equal keys are chained together
	//  into work units for the next round, using
	//  the first two keys of the group
	private size_t m_workChain;

	//  set the offset rankings and create
	//  new work units for unsorted groups
	//  of equal keys

	private void setRanks (size_t from, size_t cnt)
	{
		size_t idx = 0;

		// all members of a group get the same rank

		while( idx < cnt )
			m_rank[m_keys[from+idx++].offset] = from;

		// is this a sortable group?

		if( cnt < 2 )
			return;    // final ranking was set

		// if so, add this group to work chain for next round
		// by using the first two key prefix from the group.

		m_keys[from].prefix = m_workChain;
		m_keys[from + 1].prefix = cnt;
		m_workChain = from;
	}

	//  set the sort key (prefix) from the ranking of the offsets
	//  for rounds after the initial one.

	private void keyGroup (size_t from, size_t cnt, size_t offset)
	{
		size_t off;

		while( cnt-- ) {
			off = m_keys[from].offset + offset;
			m_keys[from++].prefix = m_rank[off];
		}
	}

	//  the tri-partite qsort partitioning

	//  creates two sets of pivot valued
	//  elements from [0:leq] and [heq:size]
	//  while partitioning a segment of the m_keys

	private void partition (size_t start, size_t size)
	{
		KeyPrefix tmp, pvt;
		KeyPrefix *lo;
		size_t loguy, higuy;
		size_t leq, heq;

		while( size > 7 ) {
			// find median-of-three element to use as a pivot
			// and swap it to the beginning of the array
			// to begin the leq group of pivot equals.

			// the larger-of-three element goes to higuy
			// the smallest-of-three element goes to middle

			lo = m_keys.ptr + start;
			higuy = size - 1;
			leq = loguy = 0;

			//  move larger of lo and hi to tmp,hi

			tmp = lo[higuy];

			if( tmp.prefix < lo.prefix )
				lo[higuy] = *lo, *lo = tmp, tmp = lo[higuy];

			//  move larger of tmp,hi and mid to hi

			if( lo[size >> 1].prefix > tmp.prefix )
				lo[higuy] = lo[size >> 1], lo[size >> 1] = tmp;

			//  move larger of mid and lo to pvt,lo
			//  and the smaller into the middle

			pvt = *lo;

			if( pvt.prefix < lo[size >> 1].prefix )
				*lo = lo[size >> 1], lo[size >> 1] = pvt, pvt = *lo;

			//  start the high group of equals
			//  with a pivot valued element, or not

			if( pvt.prefix == lo[higuy].prefix )
				heq = higuy;
			else
				heq = size;

			for( ; ; ) {
				//  both higuy and loguy are already in position
				//  loguy leaves .le. elements beneath it
				//  and swaps equal to pvt elements to leq

				while( ++loguy < higuy )
					if( pvt.prefix < lo[loguy].prefix )
						break;
					else if( pvt.prefix == lo[loguy].prefix )
						if( ++leq < loguy )
							tmp = lo[loguy], lo[loguy] = lo[leq], lo[leq] = tmp;

				//  higuy leaves .ge. elements above it
				//  and swaps equal to pvt elements to heq

				while( --higuy > loguy )
					if( pvt.prefix > lo[higuy].prefix )
						break;
					else if( pvt.prefix == lo[higuy].prefix )
						if( --heq > higuy )
							tmp = lo[higuy], lo[higuy] = lo[heq], lo[heq] = tmp;

				// quit when they finally meet at the empty middle

				if( higuy <= loguy )
					break;

				// element loguy is .gt. element higuy
				// swap them around (the pivot)

				tmp = lo[higuy];
				lo[higuy] = lo[loguy];
				lo[loguy] = tmp;
			}

			// initialize an empty pivot value group

			higuy = loguy;

			//  swap the group of pivot equals into the middle from
			//  the leq and heq sets. Include original pivot in
			//  the leq set.  higuy will be the lowest pivot
			//  element; loguy will be one past the highest.

			//  the heq set might be empty or completely full.

			if( loguy < heq )
				while( heq < size )
					tmp = lo[loguy], lo[loguy++] = lo[heq], lo[heq++] = tmp;
			else
				loguy = size;  // no high elements, they're all pvt valued

			//  the leq set always has the original pivot, but might
			//  also be completely full of pivot valued elements.

			if( higuy > ++leq )
				while( leq )
					tmp = lo[--higuy], lo[higuy] = lo[--leq], lo[leq] = tmp;
			else
				higuy = 0;    // no low elements, they're all pvt valued

			//  The partitioning around pvt is done.
			//  ranges [0:higuy-1] .lt. pivot and [loguy:size-1] .gt. pivot

			//  set the new group rank of the middle range [higuy:loguy-1]
			//  (the .lt. and .gt. ranges get set during their selection sorts)

			setRanks (start + higuy, loguy - higuy);

			//  pick the smaller group to partition first,
			//  then loop with larger group.

			if( higuy < size - loguy ) {
				partition (start, higuy);
				size -= loguy;
				start += loguy;
			} else {
				partition (start + loguy, size - loguy);
				size = higuy;
			}
		}

		//  do a selection sort for small sets by
		//  repeately selecting the smallest key to
		//  start, and pulling any group together
		//  for it at leq

		while( size ) {
			for( leq = loguy = 0; ++loguy < size; )
				if( m_keys[start].prefix > m_keys[start + loguy].prefix )
					tmp = m_keys[start], m_keys[start] = m_keys[start + loguy], m_keys[start + loguy] = tmp, leq = 0;
				else if( m_keys[start].prefix == m_keys[start + loguy].prefix )
					if( ++leq < loguy )
						tmp = m_keys[start + leq], m_keys[start + leq] = m_keys[start + loguy], m_keys[start + loguy] = tmp;

			//  now set the rank for the group of size >= 1

			setRanks (start, ++leq);
			start += leq;
			size -= leq;
		}
	}

	// the main entry point

	KeyPrefix[] sort (ubyte *buff, size_t size)
	{
		size_t start, cnt, chain;
		size_t offset = 0, off;
		size_t prefix[1];

		m_keys = new KeyPrefix[size + 1];
		m_rank = new size_t[size + prefix.sizeof];

		//  the Key and Rank arrays include stopper elements

		//m_keys = malloc ((size + 1 ) * KeyPrefix.sizeof);
		//m_keys = new KeyPrefix[size + 1].ptr;
		//memset (prefix, 0xff, prefix.sizeof);
		prefix[0] = prefix[0].max;

		// construct the suffix sorting key for each offset

		for( off = size; off--; ) {
			*prefix >>= 8;
			*prefix |= buff[off] << prefix.sizeof * 8 - 8;
			m_keys[off].prefix = *prefix;
			m_keys[off].offset = off;
		}

		// the ranking of each suffix offset,
		// plus extra ranks for the stopper elements

		//m_rank = malloc ((size + prefix.sizeof) * size_t.sizeof);

		// fill in the extra stopper ranks

		for( off = 0; off < prefix.sizeof; off++ )
			m_rank[size + off] = size + off;

		// perform the initial qsort based on the key prefix constructed
		// above.  Inialize the work unit chain terminator.

		m_workChain = size;
		partition (0, size);

		// the first pass used prefix keys constructed above,
		// subsequent passes use the offset rankings as keys

		offset = prefix.sizeof; 

		// continue doubling the key offset until there are no
		// undifferentiated suffix groups created during a run

		while( m_workChain < size ) {
			chain = m_workChain;
			m_workChain = size;

			// consume the work units created last round
			// and preparing new work units for next pass
			// (work is created in bwtsetranks)

			do {
				start = chain;
				chain = m_keys[start].prefix;
				cnt = m_keys[start + 1].prefix;
				keyGroup (start, cnt, offset);
				partition (start, cnt);
			} while( chain < size );

			//  each pass doubles the range of suffix considered,
			//  achieving Order(n * log(n)) comparisons

			offset <<= 1;
		}

		//  return the rank of offset zero in the first key

		m_keys[0].prefix = m_rank[0];
		return m_keys;
	}
}

